Activity: Area Data IV
Practice questions
Answer the following questions:
- How are row-standardized and binary spatial weights interpreted?
- What is the reason for using a Bonferroni correction for multiple tests?
- What types of spatial patterns can the local version of Moran’s I detect?
- What types of spatial patterns can the \(G_i(d)\) statistic detect?
- What is the utility of detecting hot and cold spatial spots?
Learning objectives
In this activity, you will:
- Calculate Moran’s I coefficient of autocorrelation for area data.
- Create Moran’s scatterplots.
- Examine the results of the tests/scatterplots for further insights.
- Think about ways to decide whether a landscape is random when working with area data.
Suggested reading
O’Sullivan D and Unwin D (2010) Geographic Information Analysis, 2nd Edition, Chapter 7. John Wiley & Sons: New Jersey.
Preliminaries
For this activity you will need the following:
It is good practice to clear the working space to make sure that you do not have extraneous items there when you begin your work. The command in R to clear the workspace is rm (for “remove”), followed by a list of items to be removed. To clear the workspace from all objects, do the following:
rm(list = ls())
Note that ls() lists all objects currently on the worspace.
Load the libraries you will use in this activity.
In addition to tidyverse, you will need sf, a package that implements simple features in R (you can learn about sf here) and spdep, a package that implements several spatial statistical methods (you can learn more about it here):
library(tidyverse)
library(sf)
library(spdep)
library(geog4ga3)
Begin by loading the data that you will use in this activity:
data(Hamilton_CT)
This is a sf object with census tracts and selected demographic variables for the Hamilton CMA in Canada. You can obtain new (calculated) variables as follows. For instance, to obtain the proportion of residents who are between 20 and 34 years old, and between 35 and 49:
Hamilton_CT <- mutate(Hamilton_CT, Prop20to34 = (AGE_20_TO_24 + AGE_25_TO_29 + AGE_30_TO_34)/POPULATION, Prop35to49 = (AGE_35_TO_39 + AGE_40_TO_44 + AGE_45_TO_49)/POPULATION)
You can also convert the sf object into a SpatialPolygonsDataFrame object for use with the spdedp package:
Hamilton_CT.sp <- as(Hamilton_CT, "Spatial")
This function is used to create local Moran maps:
localmoran.map <- function(p = p, listw = listw, VAR = VAR, by = by){
require(tidyverse)
require(spdep)
require(plotly)
df_msc <- transmute(p,
key = p[[by]],
Z = (p[[VAR]] - mean(p[[VAR]])) / var(p[[VAR]]),
SMA = lag.listw(listw, Z),
Type = factor(ifelse(Z < 0 & SMA < 0, "LL",
ifelse(Z > 0 & SMA > 0, "HH", "HL/LH"))))
local_I <- localmoran(p[[VAR]], listw)
df_msc <- left_join(df_msc,
data.frame(key = p[[by]], local_I))
df_msc <- rename(df_msc, p.val = Pr.z...0.)
plot_ly(df_msc) %>%
add_sf(split = ~(p.val < 0.05), color = ~Type, colors = c("red", "khaki1", "dodgerblue", "dodgerblue4"))
}
This function is used to create \(G_i^*\) maps:
gistar.map <- function(p = p, listw = listw, VAR = VAR, by = by){
require(tidyverse)
require(spdep)
require(sf)
require(plotly)
p <- mutate(p, key = p[[by]])
df.lg <- localG(p[[VAR]], listw)
df.lg <- as.numeric(df.lg)
df.lg <- data.frame(Gstar = df.lg, p.val = 2 * pnorm(abs(df.lg), lower.tail = FALSE))
df.lg <- mutate(df.lg,
Type = factor(ifelse(Gstar < 0 & p.val <= 0.05, "Low Concentration",
ifelse(Gstar > 0 & p.val <= 0.05, "High Concentration", "Not Signicant"))))
p <- left_join(p,
data.frame(key = p[[by]], df.lg))
plot_ly(p) %>%
add_sf(split = ~(p.val < 0.05), color = ~Type, colors = c("red", "dodgerblue", "gray"))
}
Create spatial weights.
- By contiguity:
Hamilton_CT.w <- nb2listw(poly2nb(pl = Hamilton_CT.sp))
- Binary, by distance (3 km threshold) including self.
Hamilton_CT.3knb <- Hamilton_CT.sp %>% coordinates() %>% dnearneigh(d1 = 0, d2 = 3)
Hamilton_CT.3kw <- nb2listw(include.self(Hamilton_CT.3knb), style = "B")
You are now ready for the next activity.
Activity
- Create local Moran maps for the population and proportion of population in the age group 20-34. What is the difference between using population (absolute) and proportion of population (rate)? Is there a reason to prefer either variable in analysis? Discuss.
data(Hamilton_CT)
pop <- mutate(Hamilton_CT, X = (AGE_20_TO_24 + AGE_25_TO_29 + AGE_30_TO_34))
proppop <- mutate(Hamilton_CT, X = (AGE_20_TO_24 + AGE_25_TO_29 + AGE_30_TO_34)/POPULATION)
pop.sp <- as(Hamilton_CT, "Spatial")
proppop.sp <- as(Hamilton_CT, "Spatial")
pop.w <- nb2listw(poly2nb(pl = pop.sp))
proppop.w <- nb2listw(poly2nb(pl = proppop.sp))
localmoran.map(pop, pop.w, "X", by = "TRACT")
Joining, by = "key"
Column `key` joining character vector and factor, coercing into character vectorNo trace type specified:
Based on info supplied, a 'scatter' trace seems appropriate.
Read more about this trace type -> https://plot.ly/r/reference/#scatter
No trace type specified:
Based on info supplied, a 'scatter' trace seems appropriate.
Read more about this trace type -> https://plot.ly/r/reference/#scatter
localmoran.map(proppop, proppop.w, "X", by = "TRACT")
Joining, by = "key"
Column `key` joining character vector and factor, coercing into character vectorNo trace type specified:
Based on info supplied, a 'scatter' trace seems appropriate.
Read more about this trace type -> https://plot.ly/r/reference/#scatter
No trace type specified:
Based on info supplied, a 'scatter' trace seems appropriate.
Read more about this trace type -> https://plot.ly/r/reference/#scatter
- Use the \(G_i^*\) statitic to analyze the population and proportion of population in the age group 20-34. What is the difference between using population (absolute) and proportion of population (rate)? Is there a reason to prefer either variable in analysis? Discuss.
gistar.map(Hamilton_CT, Hamilton_CT.3kw, "POP_DENSITY", by = "TRACT")
Joining, by = "key"
Column `key` joining character vector and factor, coercing into character vectorNo trace type specified:
Based on info supplied, a 'scatter' trace seems appropriate.
Read more about this trace type -> https://plot.ly/r/reference/#scatter
No trace type specified:
Based on info supplied, a 'scatter' trace seems appropriate.
Read more about this trace type -> https://plot.ly/r/reference/#scatter
gistar.map(Hamilton_CT, Hamilton_CT.3kw, "POPULATION", by = "TRACT")
Joining, by = "key"
Column `key` joining character vector and factor, coercing into character vectorNo trace type specified:
Based on info supplied, a 'scatter' trace seems appropriate.
Read more about this trace type -> https://plot.ly/r/reference/#scatter
No trace type specified:
Based on info supplied, a 'scatter' trace seems appropriate.
Read more about this trace type -> https://plot.ly/r/reference/#scatter
rr ggplot(Hamilton_CT, aes(x = AREA, y = POPULATION)) + geom_point()

pop <- mutate(Hamilton_CT, x = (AGE_20_TO_24 + AGE_25_TO_29 + AGE_30_TO_34))
proppop <- mutate(Hamilton_CT, x = (AGE_20_TO_24 + AGE_25_TO_29 + AGE_30_TO_34)/POPULATION)
pop.sp <- as(Hamilton_CT, "Spatial")
proppop.sp <- as(Hamilton_CT, "Spatial")
pop.w <- nb2listw(poly2nb(pl = pop.sp))
proppop.w <- nb2listw(poly2nb(pl = proppop.sp))
pop.3knb <- pop.sp %>% coordinates() %>% dnearneigh(d1 = 0, d2 = 3)
pop.3kw <- nb2listw(include.self(pop.3knb), style = "B")
proppop.3knb <- proppop.sp %>% coordinates() %>% dnearneigh(d1 = 0, d2 = 3)
proppop.3kw <- nb2listw(include.self(proppop.3knb), style = "B")
df1.lg <- localG(pop$x, pop.3kw)
df2.lg <- localG(proppop$x, proppop.3kw)
summary(df1.lg)
Min. 1st Qu. Median Mean 3rd Qu. Max.
-1.7967 -0.6468 -0.1881 0.0000 0.4402 5.2785
summary(df2.lg)
Min. 1st Qu. Median Mean 3rd Qu. Max.
-3.5782 -0.5473 -0.1022 0.0000 0.3814 3.7182
#Calcualting the p-value
df1.lg <- as.numeric(df1.lg)
df1.lg <- data.frame(Gstar = df1.lg, p.val = 2 * pnorm(abs(df1.lg), lower.tail = FALSE))
df2.lg <- as.numeric(df2.lg)
df2.lg <- data.frame(Gstar = df1.lg, p.val = 2 * pnorm(abs(df2.lg), lower.tail = FALSE))
#Join to sf:
join <- Hamilton_CT
join$Gstar <- df1.lg$Gstar
join$p.val <- df1.lg$p.val
join2 <- Hamilton_CT
join2$Gstar <- df2.lg$Gstar
join2$p.val <- df2.lg$p.val
joinplot <- mutate(join,
Type = factor(ifelse(Gstar < 0 & p.val <= 0.05, "Low Concentration",
ifelse(Gstar > 0 & p.val <= 0.05, "High Concentration", "Not Significant"))))
joinplot2 <- mutate(join2,
Type = factor(ifelse(Gstar < 0 & p.val <= 0.05, "Low Concentration",
ifelse(Gstar > 0 & p.val <= 0.05, "High Concentration", "Not Significant"))))
Error in ifelse(Gstar < 0 & p.val <= 0.05, "Low Concentration", ifelse(Gstar > :
object 'Gstar' not found
Now create local Moran maps for the population and population density in the age group 20-34. What is the difference between using population (absolute) and population density (rate)?
More generally, what do you think should guide the decision of whether to analyze variables as absolute values or rates?
LS0tDQp0aXRsZTogIkFjdGl2aXR5OiBBcmVhIERhdGEgSVYiDQpvdXRwdXQ6IGh0bWxfbm90ZWJvb2sNCi0tLQ0KDQojIEFjdGl2aXR5OiBBcmVhIERhdGEgSVYNCg0KIyMgUHJhY3RpY2UgcXVlc3Rpb25zDQoNCkFuc3dlciB0aGUgZm9sbG93aW5nIHF1ZXN0aW9uczoNCg0KMS4gSG93IGFyZSByb3ctc3RhbmRhcmRpemVkIGFuZCBiaW5hcnkgc3BhdGlhbCB3ZWlnaHRzIGludGVycHJldGVkPw0KMi4gV2hhdCBpcyB0aGUgcmVhc29uIGZvciB1c2luZyBhIEJvbmZlcnJvbmkgY29ycmVjdGlvbiBmb3IgbXVsdGlwbGUgdGVzdHM/DQozLiBXaGF0IHR5cGVzIG9mIHNwYXRpYWwgcGF0dGVybnMgY2FuIHRoZSBsb2NhbCB2ZXJzaW9uIG9mIE1vcmFuJ3MgSSBkZXRlY3Q/DQo0LiBXaGF0IHR5cGVzIG9mIHNwYXRpYWwgcGF0dGVybnMgY2FuIHRoZSAkR19pKGQpJCBzdGF0aXN0aWMgZGV0ZWN0Pw0KNS4gV2hhdCBpcyB0aGUgdXRpbGl0eSBvZiBkZXRlY3RpbmcgaG90IGFuZCBjb2xkIHNwYXRpYWwgc3BvdHM/DQoNCiMjIExlYXJuaW5nIG9iamVjdGl2ZXMNCg0KSW4gdGhpcyBhY3Rpdml0eSwgeW91IHdpbGw6DQoNCjEuIENhbGN1bGF0ZSBNb3JhbidzIEkgY29lZmZpY2llbnQgb2YgYXV0b2NvcnJlbGF0aW9uIGZvciBhcmVhIGRhdGEuDQoyLiBDcmVhdGUgTW9yYW4ncyBzY2F0dGVycGxvdHMuDQoyLiBFeGFtaW5lIHRoZSByZXN1bHRzIG9mIHRoZSB0ZXN0cy9zY2F0dGVycGxvdHMgZm9yIGZ1cnRoZXIgaW5zaWdodHMuDQozLiBUaGluayBhYm91dCB3YXlzIHRvIGRlY2lkZSB3aGV0aGVyIGEgbGFuZHNjYXBlIGlzIHJhbmRvbSB3aGVuIHdvcmtpbmcgd2l0aCBhcmVhIGRhdGEuDQoNCiMjIFN1Z2dlc3RlZCByZWFkaW5nDQoNCk8nU3VsbGl2YW4gRCBhbmQgVW53aW4gRCAoMjAxMCkgR2VvZ3JhcGhpYyBJbmZvcm1hdGlvbiBBbmFseXNpcywgMm5kIEVkaXRpb24sIENoYXB0ZXIgNy4gSm9obiBXaWxleSAmIFNvbnM6IE5ldyBKZXJzZXkuIA0KDQojIyBQcmVsaW1pbmFyaWVzDQoNCkZvciB0aGlzIGFjdGl2aXR5IHlvdSB3aWxsIG5lZWQgdGhlIGZvbGxvd2luZzoNCg0KKiBBbiBSIG1hcmtkb3duIG5vdGVib29rIHZlcnNpb24gb2YgdGhpcyBkb2N1bWVudCAodGhlIHNvdXJjZSBmaWxlKS4NCg0KKiBBIHBhY2thZ2UgY2FsbGVkIGBnZW9nNGdhM2AuDQoNCkl0IGlzIGdvb2QgcHJhY3RpY2UgdG8gY2xlYXIgdGhlIHdvcmtpbmcgc3BhY2UgdG8gbWFrZSBzdXJlIHRoYXQgeW91IGRvIG5vdCBoYXZlIGV4dHJhbmVvdXMgaXRlbXMgdGhlcmUgd2hlbiB5b3UgYmVnaW4geW91ciB3b3JrLiBUaGUgY29tbWFuZCBpbiBSIHRvIGNsZWFyIHRoZSB3b3Jrc3BhY2UgaXMgYHJtYCAoZm9yICJyZW1vdmUiKSwgZm9sbG93ZWQgYnkgYSBsaXN0IG9mIGl0ZW1zIHRvIGJlIHJlbW92ZWQuIFRvIGNsZWFyIHRoZSB3b3Jrc3BhY2UgZnJvbSBfYWxsXyBvYmplY3RzLCBkbyB0aGUgZm9sbG93aW5nOg0KYGBge3J9DQpybShsaXN0ID0gbHMoKSkNCmBgYA0KDQpOb3RlIHRoYXQgYGxzKClgIGxpc3RzIGFsbCBvYmplY3RzIGN1cnJlbnRseSBvbiB0aGUgd29yc3BhY2UuDQoNCkxvYWQgdGhlIGxpYnJhcmllcyB5b3Ugd2lsbCB1c2UgaW4gdGhpcyBhY3Rpdml0eS4gDQoNCkluIGFkZGl0aW9uIHRvIGB0aWR5dmVyc2VgLCB5b3Ugd2lsbCBuZWVkIGBzZmAsIGEgcGFja2FnZSB0aGF0IGltcGxlbWVudHMgc2ltcGxlIGZlYXR1cmVzIGluIFIgKHlvdSBjYW4gbGVhcm4gYWJvdXQgYHNmYCBbaGVyZV0oaHR0cHM6Ly9jcmFuLnItcHJvamVjdC5vcmcvd2ViL3BhY2thZ2VzL3NmL3ZpZ25ldHRlcy9zZjEuaHRtbCkpIGFuZCBgc3BkZXBgLCBhIHBhY2thZ2UgdGhhdCBpbXBsZW1lbnRzIHNldmVyYWwgc3BhdGlhbCBzdGF0aXN0aWNhbCBtZXRob2RzICh5b3UgY2FuIGxlYXJuIG1vcmUgYWJvdXQgaXQgW2hlcmVdKGh0dHBzOi8vY3Jhbi5yLXByb2plY3Qub3JnL3dlYi9wYWNrYWdlcy9zcGRlcC9pbmRleC5odG1sKSk6DQpgYGB7ciBtZXNzYWdlPUZBTFNFfQ0KbGlicmFyeSh0aWR5dmVyc2UpDQpsaWJyYXJ5KHNmKQ0KbGlicmFyeShzcGRlcCkNCmxpYnJhcnkoZ2VvZzRnYTMpDQpgYGANCg0KQmVnaW4gYnkgbG9hZGluZyB0aGUgZGF0YSB0aGF0IHlvdSB3aWxsIHVzZSBpbiB0aGlzIGFjdGl2aXR5Og0KYGBge3J9DQpkYXRhKEhhbWlsdG9uX0NUKQ0KYGBgDQoNClRoaXMgaXMgYSBgc2ZgIG9iamVjdCB3aXRoIGNlbnN1cyB0cmFjdHMgYW5kIHNlbGVjdGVkIGRlbW9ncmFwaGljIHZhcmlhYmxlcyBmb3IgdGhlIEhhbWlsdG9uIENNQSBpbiBDYW5hZGEuDQpZb3UgY2FuIG9idGFpbiBuZXcgKGNhbGN1bGF0ZWQpIHZhcmlhYmxlcyBhcyBmb2xsb3dzLiBGb3IgaW5zdGFuY2UsIHRvIG9idGFpbiB0aGUgcHJvcG9ydGlvbiBvZiByZXNpZGVudHMgd2hvIGFyZSBiZXR3ZWVuIDIwIGFuZCAzNCB5ZWFycyBvbGQsIGFuZCBiZXR3ZWVuIDM1IGFuZCA0OToNCmBgYHtyfQ0KSGFtaWx0b25fQ1QgPC0gbXV0YXRlKEhhbWlsdG9uX0NULCBQcm9wMjB0bzM0ID0gKEFHRV8yMF9UT18yNCArIEFHRV8yNV9UT18yOSArIEFHRV8zMF9UT18zNCkvUE9QVUxBVElPTiwgUHJvcDM1dG80OSA9IChBR0VfMzVfVE9fMzkgKyBBR0VfNDBfVE9fNDQgKyBBR0VfNDVfVE9fNDkpL1BPUFVMQVRJT04pDQpgYGANCg0KWW91IGNhbiBhbHNvIGNvbnZlcnQgdGhlIGBzZmAgb2JqZWN0IGludG8gYSBgU3BhdGlhbFBvbHlnb25zRGF0YUZyYW1lYCBvYmplY3QgZm9yIHVzZSB3aXRoIHRoZSBgc3BkZWRwYCBwYWNrYWdlOg0KYGBge3J9DQpIYW1pbHRvbl9DVC5zcCA8LSBhcyhIYW1pbHRvbl9DVCwgIlNwYXRpYWwiKQ0KYGBgDQoNClRoaXMgZnVuY3Rpb24gaXMgdXNlZCB0byBjcmVhdGUgbG9jYWwgTW9yYW4gbWFwczoNCmBgYHtyfQ0KbG9jYWxtb3Jhbi5tYXAgPC0gZnVuY3Rpb24ocCA9IHAsIGxpc3R3ID0gbGlzdHcsIFZBUiA9IFZBUiwgYnkgPSBieSl7DQogIHJlcXVpcmUodGlkeXZlcnNlKQ0KICByZXF1aXJlKHNwZGVwKQ0KICByZXF1aXJlKHBsb3RseSkNCiAgDQogIGRmX21zYyA8LSB0cmFuc211dGUocCwNCiAgICAgICAgICAgICAgICAgICAgICBrZXkgPSBwW1tieV1dLA0KICAgICAgICAgICAgICAgICAgICAgIFogPSAocFtbVkFSXV0gLSBtZWFuKHBbW1ZBUl1dKSkgLyB2YXIocFtbVkFSXV0pLA0KICAgICAgICAgICAgICAgICAgICAgIFNNQSA9IGxhZy5saXN0dyhsaXN0dywgWiksDQogICAgICAgICAgICAgICAgICAgICAgVHlwZSA9IGZhY3RvcihpZmVsc2UoWiA8IDAgJiBTTUEgPCAwLCAiTEwiLA0KICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIGlmZWxzZShaID4gMCAmIFNNQSA+IDAsICJISCIsICJITC9MSCIpKSkpDQogIA0KICBsb2NhbF9JIDwtIGxvY2FsbW9yYW4ocFtbVkFSXV0sIGxpc3R3KQ0KICANCiAgZGZfbXNjIDwtIGxlZnRfam9pbihkZl9tc2MsIA0KICAgICAgICAgICAgICAgICAgZGF0YS5mcmFtZShrZXkgPSBwW1tieV1dLCBsb2NhbF9JKSkNCiAgZGZfbXNjIDwtIHJlbmFtZShkZl9tc2MsIHAudmFsID0gUHIuei4uLjAuKQ0KICANCiAgcGxvdF9seShkZl9tc2MpICU+JQ0KICAgIGFkZF9zZihzcGxpdCA9IH4ocC52YWwgPCAwLjA1KSwgY29sb3IgPSB+VHlwZSwgY29sb3JzID0gYygicmVkIiwgImtoYWtpMSIsICJkb2RnZXJibHVlIiwgImRvZGdlcmJsdWU0IikpIA0KfQ0KYGBgDQoNClRoaXMgZnVuY3Rpb24gaXMgdXNlZCB0byBjcmVhdGUgJEdfaV4qJCBtYXBzOg0KYGBge3J9DQpnaXN0YXIubWFwIDwtIGZ1bmN0aW9uKHAgPSBwLCBsaXN0dyA9IGxpc3R3LCBWQVIgPSBWQVIsIGJ5ID0gYnkpew0KICByZXF1aXJlKHRpZHl2ZXJzZSkNCiAgcmVxdWlyZShzcGRlcCkNCiAgcmVxdWlyZShzZikNCiAgcmVxdWlyZShwbG90bHkpDQogIA0KICBwIDwtIG11dGF0ZShwLCBrZXkgPSBwW1tieV1dKQ0KICANCiAgZGYubGcgPC0gbG9jYWxHKHBbW1ZBUl1dLCBsaXN0dykNCiAgZGYubGcgPC0gYXMubnVtZXJpYyhkZi5sZykNCiAgZGYubGcgPC0gZGF0YS5mcmFtZShHc3RhciA9IGRmLmxnLCBwLnZhbCA9IDIgKiBwbm9ybShhYnMoZGYubGcpLCBsb3dlci50YWlsID0gRkFMU0UpKQ0KICANCiAgZGYubGcgPC0gbXV0YXRlKGRmLmxnLCANCiAgICAgICAgICAgICAgVHlwZSA9IGZhY3RvcihpZmVsc2UoR3N0YXIgPCAwICYgcC52YWwgPD0gMC4wNSwgIkxvdyBDb25jZW50cmF0aW9uIiwNCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgaWZlbHNlKEdzdGFyID4gMCAmIHAudmFsIDw9IDAuMDUsICJIaWdoIENvbmNlbnRyYXRpb24iLCAiTm90IFNpZ25pY2FudCIpKSkpDQoNCiAgcCA8LSBsZWZ0X2pvaW4ocCwgDQogICAgICAgICAgICAgICAgICBkYXRhLmZyYW1lKGtleSA9IHBbW2J5XV0sIGRmLmxnKSkNCiAgDQogIHBsb3RfbHkocCkgJT4lDQogICAgYWRkX3NmKHNwbGl0ID0gfihwLnZhbCA8IDAuMDUpLCBjb2xvciA9IH5UeXBlLCBjb2xvcnMgPSBjKCJyZWQiLCAiZG9kZ2VyYmx1ZSIsICJncmF5IikpDQp9DQpgYGANCg0KQ3JlYXRlIHNwYXRpYWwgd2VpZ2h0cy4NCg0KMSkgQnkgY29udGlndWl0eToNCmBgYHtyfQ0KSGFtaWx0b25fQ1QudyA8LSBuYjJsaXN0dyhwb2x5Mm5iKHBsID0gSGFtaWx0b25fQ1Quc3ApKQ0KYGBgDQoNCjIpIEJpbmFyeSwgYnkgZGlzdGFuY2UgKDMga20gdGhyZXNob2xkKSBfaW5jbHVkaW5nIHNlbGZfLg0KYGBge3J9DQpIYW1pbHRvbl9DVC4za25iIDwtIEhhbWlsdG9uX0NULnNwICU+JSBjb29yZGluYXRlcygpICU+JSBkbmVhcm5laWdoKGQxID0gMCwgZDIgPSAzKQ0KSGFtaWx0b25fQ1QuM2t3IDwtIG5iMmxpc3R3KGluY2x1ZGUuc2VsZihIYW1pbHRvbl9DVC4za25iKSwgc3R5bGUgPSAiQiIpDQpgYGANCg0KWW91IGFyZSBub3cgcmVhZHkgZm9yIHRoZSBuZXh0IGFjdGl2aXR5Lg0KDQojIyBBY3Rpdml0eQ0KDQoxLiBDcmVhdGUgbG9jYWwgTW9yYW4gbWFwcyBmb3IgdGhlIHBvcHVsYXRpb24gX2FuZF8gcHJvcG9ydGlvbiBvZiBwb3B1bGF0aW9uIGluIHRoZSBhZ2UgZ3JvdXAgMjAtMzQuIFdoYXQgaXMgdGhlIGRpZmZlcmVuY2UgYmV0d2VlbiB1c2luZyBwb3B1bGF0aW9uIChhYnNvbHV0ZSkgYW5kIHByb3BvcnRpb24gb2YgcG9wdWxhdGlvbiAocmF0ZSk/IElzIHRoZXJlIGEgcmVhc29uIHRvIHByZWZlciBlaXRoZXIgdmFyaWFibGUgaW4gYW5hbHlzaXM/IERpc2N1c3MuDQoNCmBgYHtyfQ0KZGF0YShIYW1pbHRvbl9DVCkNCnBvcCA8LSBtdXRhdGUoSGFtaWx0b25fQ1QsIFggPSAoQUdFXzIwX1RPXzI0ICsgQUdFXzI1X1RPXzI5ICsgQUdFXzMwX1RPXzM0KSkNCnByb3Bwb3AgPC0gbXV0YXRlKEhhbWlsdG9uX0NULCBYID0gKEFHRV8yMF9UT18yNCArIEFHRV8yNV9UT18yOSArIEFHRV8zMF9UT18zNCkvUE9QVUxBVElPTikNCnBvcC5zcCA8LSBhcyhIYW1pbHRvbl9DVCwgIlNwYXRpYWwiKQ0KcHJvcHBvcC5zcCA8LSBhcyhIYW1pbHRvbl9DVCwgIlNwYXRpYWwiKQ0KcG9wLncgPC0gbmIybGlzdHcocG9seTJuYihwbCA9IHBvcC5zcCkpDQpwcm9wcG9wLncgPC0gbmIybGlzdHcocG9seTJuYihwbCA9IHByb3Bwb3Auc3ApKQ0KDQpsb2NhbG1vcmFuLm1hcChwb3AsIHBvcC53LCAiWCIsIGJ5ID0gIlRSQUNUIikNCmxvY2FsbW9yYW4ubWFwKHByb3Bwb3AsIHByb3Bwb3AudywgIlgiLCBieSA9ICJUUkFDVCIpDQpgYGANCg0KMi4gVXNlIHRoZSAkR19pXiokIHN0YXRpdGljIHRvIGFuYWx5emUgdGhlIHBvcHVsYXRpb24gX2FuZF8gcHJvcG9ydGlvbiBvZiBwb3B1bGF0aW9uIGluIHRoZSBhZ2UgZ3JvdXAgMjAtMzQuIFdoYXQgaXMgdGhlIGRpZmZlcmVuY2UgYmV0d2VlbiB1c2luZyBwb3B1bGF0aW9uIChhYnNvbHV0ZSkgYW5kIHByb3BvcnRpb24gb2YgcG9wdWxhdGlvbiAocmF0ZSk/IElzIHRoZXJlIGEgcmVhc29uIHRvIHByZWZlciBlaXRoZXIgdmFyaWFibGUgaW4gYW5hbHlzaXM/IERpc2N1c3MuDQoNCmBgYHtyfQ0KZ2lzdGFyLm1hcChIYW1pbHRvbl9DVCwgSGFtaWx0b25fQ1QuM2t3LCAiUE9QX0RFTlNJVFkiLCBieSA9ICJUUkFDVCIpDQpnaXN0YXIubWFwKEhhbWlsdG9uX0NULCBIYW1pbHRvbl9DVC4za3csICJQT1BVTEFUSU9OIiwgYnkgPSAiVFJBQ1QiKQ0KYGBgDQoNCmBgYHtyfQ0KZ2dwbG90KEhhbWlsdG9uX0NULCBhZXMoeCA9IEFSRUEsIHkgPSBQT1BVTEFUSU9OKSkgKyBnZW9tX3BvaW50KCkNCmBgYA0KDQoNCg0KYGBge3J9DQpwb3AgPC0gbXV0YXRlKEhhbWlsdG9uX0NULCB4ID0gKEFHRV8yMF9UT18yNCArIEFHRV8yNV9UT18yOSArIEFHRV8zMF9UT18zNCkpDQpwcm9wcG9wIDwtIG11dGF0ZShIYW1pbHRvbl9DVCwgeCA9IChBR0VfMjBfVE9fMjQgKyBBR0VfMjVfVE9fMjkgKyBBR0VfMzBfVE9fMzQpL1BPUFVMQVRJT04pDQpwb3Auc3AgPC0gYXMoSGFtaWx0b25fQ1QsICJTcGF0aWFsIikNCnByb3Bwb3Auc3AgPC0gYXMoSGFtaWx0b25fQ1QsICJTcGF0aWFsIikNCnBvcC53IDwtIG5iMmxpc3R3KHBvbHkybmIocGwgPSBwb3Auc3ApKQ0KcHJvcHBvcC53IDwtIG5iMmxpc3R3KHBvbHkybmIocGwgPSBwcm9wcG9wLnNwKSkNCg0KcG9wLjNrbmIgPC0gcG9wLnNwICU+JSBjb29yZGluYXRlcygpICU+JSBkbmVhcm5laWdoKGQxID0gMCwgZDIgPSAzKQ0KcG9wLjNrdyA8LSBuYjJsaXN0dyhpbmNsdWRlLnNlbGYocG9wLjNrbmIpLCBzdHlsZSA9ICJCIikNCg0KcHJvcHBvcC4za25iIDwtIHByb3Bwb3Auc3AgJT4lIGNvb3JkaW5hdGVzKCkgJT4lIGRuZWFybmVpZ2goZDEgPSAwLCBkMiA9IDMpDQpwcm9wcG9wLjNrdyA8LSBuYjJsaXN0dyhpbmNsdWRlLnNlbGYocHJvcHBvcC4za25iKSwgc3R5bGUgPSAiQiIpDQoNCmRmMS5sZyA8LSBsb2NhbEcocG9wJHgsIHBvcC4za3cpDQpkZjIubGcgPC0gbG9jYWxHKHByb3Bwb3AkeCwgcHJvcHBvcC4za3cpDQpzdW1tYXJ5KGRmMS5sZykNCnN1bW1hcnkoZGYyLmxnKQ0KDQojQ2FsY3VhbHRpbmcgdGhlIHAtdmFsdWUgDQpkZjEubGcgPC0gYXMubnVtZXJpYyhkZjEubGcpDQpkZjEubGcgPC0gZGF0YS5mcmFtZShHc3RhciA9IGRmMS5sZywgcC52YWwgPSAyICogcG5vcm0oYWJzKGRmMS5sZyksIGxvd2VyLnRhaWwgPSBGQUxTRSkpDQoNCmRmMi5sZyA8LSBhcy5udW1lcmljKGRmMi5sZykNCmRmMi5sZyA8LSBkYXRhLmZyYW1lKEdzdGFyID0gZGYxLmxnLCBwLnZhbCA9IDIgKiBwbm9ybShhYnMoZGYyLmxnKSwgbG93ZXIudGFpbCA9IEZBTFNFKSkNCg0KI0pvaW4gdG8gc2Y6IA0Kam9pbiA8LSBIYW1pbHRvbl9DVA0Kam9pbiRHc3RhciA8LSBkZjEubGckR3N0YXINCmpvaW4kcC52YWwgPC0gZGYxLmxnJHAudmFsIA0Kam9pbjIgPC0gSGFtaWx0b25fQ1QNCmpvaW4yJEdzdGFyIDwtIGRmMi5sZyRHc3Rhcg0Kam9pbjIkcC52YWwgPC0gZGYyLmxnJHAudmFsIA0KDQpqb2lucGxvdCA8LSBtdXRhdGUoam9pbiwgDQogICAgICAgICAgICAgICAgICAgVHlwZSA9IGZhY3RvcihpZmVsc2UoR3N0YXIgPCAwICYgcC52YWwgPD0gMC4wNSwgIkxvdyBDb25jZW50cmF0aW9uIiwNCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBpZmVsc2UoR3N0YXIgPiAwICYgcC52YWwgPD0gMC4wNSwgIkhpZ2ggQ29uY2VudHJhdGlvbiIsICJOb3QgU2lnbmlmaWNhbnQiKSkpKQ0Kam9pbnBsb3QyIDwtIG11dGF0ZShqb2luMiwgDQogICAgICAgICAgICAgICAgICAgVHlwZSA9IGZhY3RvcihpZmVsc2UoR3N0YXIgPCAwICYgcC52YWwgPD0gMC4wNSwgIkxvdyBDb25jZW50cmF0aW9uIiwNCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBpZmVsc2UoR3N0YXIgPiAwICYgcC52YWwgPD0gMC4wNSwgIkhpZ2ggQ29uY2VudHJhdGlvbiIsICJOb3QgU2lnbmlmaWNhbnQiKSkpKQ0KDQojcGxvdA0KZ2dwbG90KGRhdGEgPSBqb2lucGxvdCkgKw0KICBnZW9tX3NmKGFlcyhmaWxsID0gVHlwZSkpDQpnZ3Bsb3QoZGF0YSA9IGpvaW5wbG90MikgKw0KICBnZW9tX3NmKGFlcyhmaWxsPVR5cGUpKQ0KDQpgYGANCg0KDQoNCjMuIE5vdyBjcmVhdGUgbG9jYWwgTW9yYW4gbWFwcyBmb3IgdGhlIHBvcHVsYXRpb24gX2FuZF8gcG9wdWxhdGlvbiBkZW5zaXR5IGluIHRoZSBhZ2UgZ3JvdXAgMjAtMzQuIFdoYXQgaXMgdGhlIGRpZmZlcmVuY2UgYmV0d2VlbiB1c2luZyBwb3B1bGF0aW9uIChhYnNvbHV0ZSkgYW5kIHBvcHVsYXRpb24gZGVuc2l0eSAocmF0ZSk/DQoNCjQuIE1vcmUgZ2VuZXJhbGx5LCB3aGF0IGRvIHlvdSB0aGluayBzaG91bGQgZ3VpZGUgdGhlIGRlY2lzaW9uIG9mIHdoZXRoZXIgdG8gYW5hbHl6ZSB2YXJpYWJsZXMgYXMgYWJzb2x1dGUgdmFsdWVzIG9yIHJhdGVzPw==